NASA/TM-2009-2 15705 



Compressibility Considerations for k-cq 
Turbulence Models in Hypersonic Boundary 
Layer Applications 


C. L. Rumsey 

Langley Research Center, Hampton, Virginia 


April 2009 



NASA STI Program ... in Profile 


Since its founding, NASA has been dedicated to 
the advancement of aeronautics and space science. 

The NASA scientific and technical information (STI) 
program plays a key part in helping NASA maintain 
this important role. 

The NASA STI program operates under the 
auspices of the Agency Chief Information Officer. It 
collects, organizes, provides for archiving, and 
disseminates NASA’s STI. The NASA STI program 
provides access to the NASA Aeronautics and Space 
Database and its public interface, the NASA Technical 
Report Server, thus providing one of the largest 
collections of aeronautical and space science STI in 
the world. Results are published in both non-NASA 
channels and by NASA in the NASA STI Report 
Series, which includes the following report types: 

• TECHNICAL PUBLICATION. Reports of 
completed research or a major significant phase 
of research that present the results of NASA 
programs and include extensive data or 
theoretical analysis. Includes compilations of 
significant scientific and technical data and 
infonnation deemed to be of continuing 
reference value. NASA counterpart of peer- 
reviewed formal professional papers, but having 
less stringent limitations on manuscript length 
and extent of graphic presentations. 

• TECHNICAL MEMORANDUM. Scientific 
and technical findings that are preliminary or of 
specialized interest, e.g., quick release reports, 
working papers, and bibliographies that contain 
minimal annotation. Does not contain extensive 
analysis. 

• CONTRACTOR REPORT. Scientific and 
technical findings by NASA-sponsored 
contractors and grantees. 


• CONFERENCE PUBLICATION. Collected 
papers from scientific and technical 
conferences, symposia, seminars, or other 
meetings sponsored or co-sponsored by NASA. 

• SPECIAL PUBLICATION. Scientific, 
technical, or historical information from NASA 
programs, projects, and missions, often 
concerned with subjects having substantial 
public interest. 

• TECHNICAL TRANSLATION. English- 
language translations of foreign scientific and 
technical material pertinent to NASA’s mission. 

Specialized services also include creating custom 

thesauri, building customized databases, and 

organizing and publishing research results. 

For more information about the NASA STI 

program, see the following: 

• Access the NASA STI program home page at 
http://www.sti. nasa.gov 

• E-mail your question via the Internet to 
help@sti.nasa.gov 

• Fax your question to the NASA STI Help Desk 
at 443-757-5803 

• Phone the NASA STI Help Desk at 
443-757-5802 

• Write to: 

NASA STI Help Desk 

NASA Center for AeroSpace Information 

7115 Standard Drive 

Hanover, MD 21076-1320 


NASA/TM-2009-2 15705 



Compressibility Considerations for k-cq 
Turbulence Models in Hypersonic Boundary 
Layer Applications 


C. L. Rumsey 

Langley Research Center, Hampton, Virginia 


National Aeronautics and 
Space Administration 

Langley Research Center 
Hampton, Virginia 23681-2199 


April 2009 



Available from: 

NASA Center for AeroSpace Information 
7115 Standard Drive 
Hanover, MD 21076-1320 
443-757-5802 





Abstract 


The ability of k-u models to predict compressible turbulent skin friction in hypersonic boundary 
layers is investigated. Although uncorrected two-equation models can agree well with correlations 
for hot-wall cases, they tend to perform progressively worse - particularly for cold walls - as the 
Mach number is increased in the hypersonic regime. Simple algebraic models such as Baldwin- 
Lomax perform better compared to experiments and correlations in these circumstances. Many of 
the compressibility corrections described in the literature are summarized here. These include cor- 
rections that have only a small influence for k-u models, or that apply only in specific circumstances. 
The most widely-used general corrections were designed for use with jet or mixing-layer free shear 
flows. A less well-known dilatation-dissipation correction intended for boundary layer flows is also 
tested, and is shown to agree reasonably well with the Baldwin-Lomax model at cold-wall condi- 
tions. It exhibits a less dramatic influence than the free shear type of correction. There is clearly a 
need for improved understanding and better overall physical modeling for turbulence models applied 
to hypersonic boundary layer flows. 


1 Introduction 


Compressibility is typically not considered to be important for wall-bounded turbulent flows over a 
wide range of Mach numbers. As stated in Wilcox [1] (p. 239): “Generally speaking, compressibil- 
ity has a relatively small effect on turbulent eddies in wall-bounded flows. This appears to be true 
for Mach numbers up to about 5 (and perhaps as high as 8), provided the flow doesn’t experience 
large pressure changes over a short distance such as we might have across a shock wave. At sub- 
sonic speeds, compressibility effects on eddies are usually unimportant for boundary layers provided 
T w /T e < 6.” 

The hypothesis of Morkovin [2] states that the compressibility effects on turbulence can be ac- 
counted for by mean density variations alone. For many applications, this hypothesis has proved 
correct in that good results can be obtained for mean velocity and temperature fields using incom- 
pressible turbulence models extended directly to compressible turbulent boundary layers. Further- 
more, So et al. [3] have shown the Morkovin hypothesis to be equally applicable for prediction of 
the turbulence field itself, for flat plate boundary layers up to a Mach number of at least 10. They 
state: “there is indeed a dynamic similarity of the incompressible and compressible mean and turbu- 
lence field, and the Morkovin hypothesis is valid for both fields.” In other words, for many subsonic 
through hypersonic boundary layer applications, the incompressible forms of turbulence models 
(with mean density variations accounted for) are expected to be reasonable approximations. 

The most common classes of compressibility correction for Reynolds-averaged Navier-Stokes 
(RANS) turbulence models were developed for the purpose of improving correlations with experi- 
ment for free shear layer or jet spreading rates. See, e.g.. Refs. 4-6. However, what we are concerned 
with here is primarily (attached) hypersonic boundary layer flow. In this paper, compressibility cor- 
rections (particularly applicable to boundary layer flows) from the literature are described. The 
focus here is solely on the k-u form of two-equation models. The claim that compressibility correc- 
tions are not required for hypersonic boundary layer flows is investigated for a wide range of Mach 
numbers and wall-temperature boundary conditions. 

The paper is organized as follows. First, several standard forms of the k-u model are given. 
Then, compressibility corrections from the literature are described. Finally, results for hypersonic 
boundary layer flows are shown, and conclusions are made. 
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2 Standard Forms of Two-Equation Turbulence Models 

2.1 Wilcox 1988 

The original incompressible form of the Wilcox k-to model [7], referred to here as Wilcox88, is 
written as: 
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Note that the definition for r- i? - varies in the literature: sometimes it is defined with the opposite 
sign, and sometimes it is defined without the density. The definition does not matter as long as the 
production term is defined appropriately in Eq. (4), with +2fi t Sij(dui/dxj ) as the leading term in 
V. In Eq. (5), the — (l/3)(du k /dx k )8ij term and the (2/3 )pk5ij term are often ignored for low- 
speed flows (the former term makes the strain rate tensor traceless in 3-D flows), but these both may 
be non-negligible for higher speed flows, or near stagnation regions. The constants are (3* = 0.09, 
a k = 0.5, 7 = 5/9, /3 = 3/40, and a u = 0.5. 

For clarity, the production term is expanded out here: 
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where ,S', , is the traceless form of the strain rate tensor (in 3-D): 
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2.2 Wilcox 2006 

In 1998 Wilcox presented a modified version of the k-u> model [ 8 ]. Because it has been superseded, 
this 1998 version is not described here. A newer form of the Wilcox k-u model [1,9], referred to 
here as Wilcox06, was developed to improve the predictive accuracy compared to the 1988 and 1998 
versions for free shear flows and strongly separated flows (and hence be more competitive with the 
Menter SST model, described in the next section). Wilcox06 is given by: 
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and Cum = 7/8. The coefficients are 7 = 13/25, (3 = (3$ fp, (3* = 0.09, a w = 0.5, a k = 3/5, 
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where Ski = Ski — \{du m / dx m )5ki- This form forces Xu> = 0 for 2-D flow (both incompressible 
and compressible). The f rS parameter was added by Wilcox to account for the round-jet/plane-jet 
anomaly [1], For boundary layer applications, the vortex stretching parameter is sometimes 
ignored (set to zero), yielding fp= 1. 


2.3 Menter SST 

The two-equation SST model of Menter [10] is written as: 
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The eddy viscosity is given by: 
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where a± = 0.31, U is the magnitude of vorticity, and F\ and f 2 are blending functions (given 
below). 

The “shear stress transport” (SST) part of the model is based on Bradshaw’s assumption that 
the principal shear stress is proportional to k, via: T12 = —pa\k. From Eq. (5), the primary term 
in eddy viscosity models is: T12 = — 2 p t Si 2 - In adverse pressure gradient boundary layer flows, 
the standard method often leads to too much eddy viscosity (an overprediction of T12), inhibiting or 
delaying separation. In these situations, it is better for the model to choose T12 based on Bradshaw’ 
assumption. Using: 


paik = 2p t Si2 (25) 

we find how to set the eddy viscosity in order to recover values corresponding with Bradshaw’s 
assumption: 
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Functionally, Eq. (24) chooses the minimum eddy viscosity between the standard one and that dic- 
tated by Bradshaw’s assumption limited to within the boundary layer region. 

In the SST model, there are two sets of coefficients, which are combined using a blending 
function. The constants for set 1 are (3\ = 0.09, <j k 1 = 0.85, (3i = 0.075, <t w i = 0.5, and 
7i = di/Pi ~ <Xljik 2 / sfdi ~ 0.55317. The constants for set 2 are /?£ = 0.09, Ofc 2 = 1.0, 
P 2 = 0.0828, = 0.856, and 72 = P 2 /P 2 ~ ~ 0.44035. The constant k is defined 

as k = 0.41. Set 1 and set 2 are blended via: 
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where d is the distance to the nearest wall. The b\ term is given by: 


F 2 = tanh (arg 2 ) 


(31) 
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2.4 Other Considerations for k-oj Models 

For situations in which compressibility is important (and shocks may be present), the turbulence 
equations are sometimes solved in conservation form. For example, making use of the continuity 
equation dp/dt + d(puj)/dxj = 0, Eqs. (1) and (2) can be written: 
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It is unclear whether the turbulence model equation form (conservative vs. non-conservative) makes 
much difference in the common situation where the turbulence models are solved separately (loosely 
coupled) from the conservative mean flow equations. But certainly if the equations are fully coupled, 
all should be solved consistently in conservation form. 

For flows in which the turbulent kinetic energy is non-negligible compared to the square of the 
mean velocity, the k contributes to the conservation of total energy via pE = p(e + \uiUi + k) (see 
Wilcox [1]). Also, the molecular and turbulent diffusion of /,:, typically modeled as 
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dxj 
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in the mean flow energy equation [11], should be included. Furthermore, the perfect gas equation of 
state then becomes [12]: p = (7 — 1 )p(E — \uiUi — k). Because many CFD codes include other 
(simpler) turbulence models beside two-equation models, for which k is not available, the turbulent 
kinetic energy contribution to total energy (and its explicit appearance in the energy equation and 
equation of state) is often ignored. 

In an often-used variant of the k-oj model, the production term is simplified by an approximation 
that makes use of the local magnitude of vorticity O: 


2 du.j q 2 du k 

V = p t ii - -pkdij - — = p t fl 2 3 - -pk- — (36) 

3 dxj 3 dx k 

This vorticity source term is often a good approximation of the exact source term in boundary layer 
flows [13], and its use can avoid some numerical difficulties sometimes associated with the use of 
the exact source term. Again, it is often common to ignore the (2/3) pk term in the production source 
of Eq. (36) for many applications, but this may have a non-negligible influence for high-speed flows. 

The recommended wall boundary conditions for k-ui models are [10]: k wa ii = 0, co wa u = 
60i'/[/3i(A(ii) 2 ]. In his boundary layer code, Wilcox [1] overwrites the computed value of ui with 
the theoretical value u> « 6 ^/ (/?i y 2 ) at the first grid point off the surface (for smooth walls), but this 
method of overwriting field variables rather than specifying boundary conditions is undesirable in 
general-purpose Navier-Stokes codes. The farfield boundary conditions are more difficult to define 
with confidence. Part of the problem is that the freestream levels are not preserved; they decay 
rapidly (both due to the equations themselves as well as due to typically coarse grid spacing in the 
farfield). This decay, which occurs for k-e equations as well, makes the local “ambient” levels near 
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the body a function of the farfield grid extent. In the freestream, the k-u> governing equations dictate 
that the decay of eddy viscosity occurs according to: 


l^t — /A,oo 


1 + /JWqo 
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1 — ( 0 . 09 / ( 3 ) 


(37) 


where x is the distance from the location where the boundary conditions are set. As discussed 
in Spalart and Rumsey [14], real flow over external aerodynamic configurations has no reason to 
obey the decay equations used to calibrate two-equation models in isotropic turbulence. In reality, 
the kinetic energy (and eddy viscosity) relevant to the aircraft flow varies very little over the size 
of the typical CFD domain. Thus, the behavior represented by decaying freestream turbulence is 
not representative of reality. Ref. 14 describes the use of sustaining terms which, when added to 
the k-oj equations, preserve the freestream levels without decay. Although maintaining freestream 
turbulence is probably important from the point of view of numerical-transition consistency with 
grid refinement, at high Reynolds numbers the effect is generally small. Therefore, results in this 
paper will not use the sustaining terms. 

For typical subsonic/transonic/supersonic applications, most CFD application codes have devel- 
oped their own methodology for setting farfield boundary conditions for k and u>, in order to yield 
reasonable results across a broad range of problems. For example, in CFL3D (Krist et al. [15]), the 
boundary conditions are: k/a ^ = 9 x 10 -9 and u}n 00 /(p 00 a^ X) ) = 1 x 10~ 6 , which always gives 
p t , oo/Voo = 0.009. Because the freestream turbulence level, Tu (in percent), is given by: 


Tu, % = 100 



(38) 


this means that a fixed k/a ^ = 9 x 10 9 yields (for example): Tu = 0.0387% for M = 0.2, 
0.0039% for M = 2.0, 0.0016% for M = 5.0, and 0.0008% for M = 10.0. Thus, perhaps, it 
makes better sense when solving over a broad range of Mach numbers to use a fixed Tu (i.e., fixed 
k/u/^) in the freestream instead of fixed k/a ^,). Otherwise higher Mach number cases will have an 
even greater tendency to become laminar. 


3 Compressibility Corrections 

Wilcox [1] describes many of the compressible-flow closure approximations. A few of them are 
already employed in most compressible flow CFD codes. For example, the Reynolds stress tensor of 
Eq. (5) is already written appropriately for compressible flows. The most commonly used turbulent 
heat-flux vector (qx = — ( p t J Pr t )dh/ dxj , where h is enthalpy and Pr t is typically around 0.9 
for boundary layers), has been in common use in compressible flow CFD codes for many years. 
However, the models for pressure-diffusion, pressure-dilatation, and pressure-work are all either 
under development, very little is known, or proposed models are too complex or have not gained 
wide acceptance (see, e.g., Zeman [16], Grasso and Falconi [17], and Yoshizawa et al. [18]). Many 
of these compressibility effects are presumed to be small in boundary layers [1]. As a result, most 
widely-used models do not include them. For example, Sarkar’s model for the pressure-dilatation 
correction in compressible flows [19] is rarely employed for boundary-layer computations. See 
also Wilcox [1] and Grasso and Falconi [17]. In the Sarkar model, the pressure -dilatation adds the 
following term to the /.'-equation (in the k-e model): 

{- 0 . 2 V + a 3 pe)Mx (39) 

where «2 = 0-15 (0.4 in Ref. 17), a 3 = 0.2, Mt = (V2k) /a, and a is the local speed of sound. 

On the other hand, the Sarkar/Zeman compressibility corrections for dilatation-dissipation are 
often employed for jets and free shear mixing layers, in spite of the fact that the reasoning behind 
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them is fundamentally flawed [1], See also Sarkar [20]. It turns out that dilatation-dissipation 
is small or negligible, and mixing-layer compressibility effects likely manifest themselves in the 
pressure-strain redistribution term. Properly formulated corrections are still being explored. In 
the mean time, the existing dilatation-dissipation corrections provide the desired trends for mixing 
layers, albeit for the wrong reasons, and so are still considered useful for those cases. It should also 
be noted that dilatation-dissipation models typically account for the pressure-dilatation correction 
[21], so when employing both corrections the coefficient in the dilatation-dissipation model must be 
reduced by about a factor of two from its standard value [17]. 

Unfortunately, the Sarkar/Zeman compressibility corrections can have a detrimental effect on 
many boundary-layer predictions (they tend to produce wall skin frictions that are too low, and can 
also negatively impact the size of the predicted separation region [22]). Wilcox [1,23] developed a 
modification that significantly decreases this detrimental effect, and Brown [24] further attempted to 
eliminate its potential impact in very high Mach number boundary layers by combining it with the 
Fi function of Menter (Eq. (28)). In the Wilcox correction, the coefficients of the k-ui destruction 
terms are modified as follows: 

/3* = p* [1 + CF(M t )} (40) 

/3c = P-P*CF{M t ) (41) 

where 

F{M t ) = (M| - M f, o ) H (M r - M Tq ) (42) 

and = 2, M To = 0.25, M T = (V2 k)/a, H(-) is the Heaviside function, and a is the local speed 
of sound. 

However, it has been observed 1 that for cold-wall cases the skin friction is typically overpre- 
dicted to such an extent that including a dilatation-dissipation correction can yield improved results, 
although possibly for the wrong reasons. Zeman [16] noted that the apparent unimportance of the 
pressure-dilatation and dilatation-dissipation in boundary layers is “only a question of degree.” He 
found that as freestream Mach number increases and wall cooling increases, compressibility effects 
become increasingly important. In the Zeman dilatation-dissipation correction for boundary layers, 
the coefficients of the k-uj destruction terms are modified as in Eqs. (40) and (41), only now F(Mt) 
is given by: 


F{M t ) 


1 — exp 


l' Mt — Mt 0 

V A 


j H (M t — M To ) 


(43) 


with = 0.75, M To = 0.2, and A = 0.66. 

In addition to a pressure-dilatation correction and a dilatation-dissipation correction, Grasso and 
Falconi [17] also included a correction to the A: -equation in their k-e model, due to the scalar product 
of the Favre-velocity and the mean pressure gradient. They believed that this term may be influential 
in regions of large pressure and density gradients. 

Both Wilcox [1] and Huang et al. [25] mention that the k-e form of the two-equation model 
exhibits more deviation from the compressible law-of-the-wall than k-oj at high Mach numbers. 
Wilcox also points out that the k-e form is more problematic for adverse pressure-gradient wall- 
bounded flows. Huang et al. proposed a possible iterative procedure to reproduce the expected 
profile, and also mentioned that alternative forms such as k-(e 5 ^ e /k) may reduce the sensitivity, but 
neither of these proposals were widely used. 

Catris et al. [26] extended the analysis of Huang et al. [25], and showed that specific corrections 
to the diffusion terms are necessary to make the models consistent with the logarithmic law for 

'White, J. A., private communication 2008. 
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compressible boundary layers. The corrections were derived for a variety of models. Here, we are 
only concerned with the k-oj form. For example, Eqs. (33) and (34) get altered as follows: 
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where the only changes are in the diffusion terms. However, Catris et al. point out that for k-u>, the 
difference in results due to modifying the diffusion terms is only very slight. This further confirms 
the low sensitivity of k-u> to compressibility effects in boundary layer flow, as described in Refs. 1 
and 25. 

A length-scale modification was proposed by Vuong and Coakley [27] and Huang and Coak- 
ley [28], to reduce the magnitude of heat transfer for high-speed separated boundary layers near 
reattachment. See also Brown [24] and Coratekin et al. [29]. In this correction, the length scale 
going into the eddy viscosity is limited based on Bradshaw’s relation r t /p oc a±k and the mixing 
length relation £ = k d. The turbulent length scale is limited by: 


£ = min 





(46) 


(Note that in Vuong and Coakley and Huang and Coakley, there is an additional factor of 6), = 0.09 
present due to the different way that u> is defined.) Because u> = \/k/£, the end result is that the 
eddy viscosity is limited according to: 


Pt 
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where pt. tS td is the eddy viscosity computed the usual way, via Eq. (3), (16), or (24). It should be 
noted here that the assumption of constant turbulent Prandtl number has been recently questioned, 
relative to its effect on heat flux for shock/boundary layer cases [30], A variable turbulent Prandtl 
number model was shown to improve heat flux near reattachment. 

A rapid compression fix was implemented by Coakley and Huang [31] (see also Vuong and 
Coakley [27], Coratekin et al. [29], and Forsythe et al. [32]). In this fix, the production term in the 
^-equation 
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or equivalently: 



where ,S/ 7 was defined earlier in Eq. (13). Thus, in this rapid compression fix, the original (2/3)7 
(which is close to 1/3) gets increased to a fixed value of 4/3 for linear deformations. 

This modification was made in order to increase the size of computed separation bubble regions, 
by insuring that the turbulent length scale does not change too quickly when undergoing rapid com- 
pression. However, as discussed in Forsythe et al. [32], the shear-stress transport part of the Menter 
SST model already improves correlations with experimental bubble size, so the ad hoc rapid com- 
pression fix was not used for that model. The Wilcox06 model is designed with a similar stress 
limiter modification, so this model, too, probably would not benefit from the rapid compression fix. 

In summary, when considering high Mach number compressible boundary layer flows using 
k-oj models, the conservation of total energy should be configured to include the contribution of 
the turbulent kinetic energy k, and the mean flow energy equation should include the molecular 
and turbulent diffusion of k. It is sometimes common practice to ignore these effects, which is 
certainly justified when k is significantly smaller than the square of the mean velocity. Furthermore, 
turbulent production should officially include the (2/3)pfc term, which multiplies duk/dxu and 
hence is identically zero for incompressible flows. Because this term typically has little effect over a 
broad range of conditions, it is sometimes ignored, particularly when other approximations, limiting, 
etc. are employed for turbulence production. Other than these considerations (which may or may 
not be important, depending on the case), little extra appears to be called for - based on the currently 
available literature - in terms of specific corrections for compressible Navier-Stokes codes applied 
to boundary layer flows. Adopting the modified diffusion-term form of Catris et al. [26] has been 
shown to make very little difference for k-u models. The length-scale modification of Vuong and 
Coakley [27] appears to only be important for the specific circumstance of predicting heat transfer 
near reattachment after separation, but this flow feature may also be improved by adopting models 
with variable turbulent Prandtl number. The rapid compression fix of Coakley and Huang [31] has 
been negated by the more accepted stress limiters that appear in SST and Wilcox06. 

However, Zeman [16] and practical experience indicates that the need for corrections in hyper- 
sonic boundary layers becomes increasingly evident as Mach number increases, particularly for cold 
walls. But the “traditional” Sarkar / Zeman / Wilcox fixes for free shear flows tend to over-correct 
in many cases when applied to boundary layer flows. In the Results section, the influence of the 
less-widely-used Zeman correction (formulated for boundary layer flows) is explored. 


4 Results 

4.1 Experimental Correlations for Skin Friction and Heat Transfer on a Flat 
Plate 

Wall skin friction is given by the formula: 


Cf = 


T~w 

2peU? 


(53) 


where t w is the wall shear stress y w du/dy\ w and the subscript “e” represents “edge” or freestream 
values. There have been a plethora of correlations for wall skin friction (and heat transfer) for the 
flat plate over the years. See, for example. White [33], Peterson [34], and Hopkins and Inouye [35], 
One of the reasons for the large number is the fact that there is a significant amount scatter in the 
available experimental data, especially those with heat transfer, making certainty difficult. Many of 
the skin friction correlations use a compressibility transformation idea: 


Cy ^ f ^incompi-^^O ^Ree') 

■F c 


( 54 ) 



In other words, the formula for incompressible Cf/ ulcoul]) , which is often expressed as a function of 
Reg, is instead computed using the altered variable RegFp> ee , and then divided by the function F c . 
A widely-used correlation for Cf t i ncom p is the Karman-Schoenherr relation (see Roy and Blottner 
[36]): 


C/, incomp logio (2i?ee) [I7.0751og 10 ( 2 Reg) + 14.832] ? ' 

Here, we examine three correlations for Cf \ van Driest [37], Spalding and Chi [38], and White and 
Christoph [33], For each of these, the F c is defined the same way: 


F r = 


- 1 - aw / e 


(sin 1 A + sin 1 B) 
where T e represents the “edge” or freestream temperature, and 


(56) 


T aw =T e (l + r M 2 ^J (57) 

The recovery factor r is taken to be 0.9. This empirical factor is often introduced because in practice 
energy recovery is not perfect. The numerator of Eq. (56) is thus simply |r ( 7 — 1 )M 2 . The A and 
B are given by: 


2 a 2 — b 
{b 2 + 4a 2 ) 1 / 2 


{b 2 + 4a 2 ) 1 / 2 


where 


(58) 

(59) 


7-1 2 T e \ 1/2 

r ~ M t) m 

* = (<>i) 

A contour plot of F c as a function of Mach number and log 10 (T w /T e ) is shown in Fig. 1 . A different 
version of this same plot is also given in Spalding and Chi [38]. 

The three correlations differ in their definitions of Fn eg : 




\ 0.772 

I Spalding and Chi 

1/2 

White and Christoph 


(62) 

(63) 

(64) 


(Note the typographical error in White [33] in Table 7-3, where the term Q is inverted.) Both van 
Driest and White/Christoph correlations are functions of n e //j, w . This quantity can be obtained via 
Sutherland’s law (see White [33]): 


M = Mo 


f T \ 

\To) 


T 0 + S' 


(65) 


10 


T + S' 



where for air /i 0 = 0.1716 mP, T 0 = 491.6 R, and S' = 199 R. Thus: 

ife = _3/2 (r w /T e ) + (syr e ) 

Mu> V^j 1 + (S'/T e ) 1 ; 

Therefore, when Sutherland’s law is used, both of these correlations are functions not only of the 
ratio T w /T e , but also of the freestream temperature T e itself. For all of the work herein, T e is chosen 
to be 540 R. Contour plots of Fn ee for the three correlations are shown in Figs. 2(a) - (c). Note that 
a different version of Fig. 2(b) can also be found in Spalding and Chi [38]. 

Using Eq. (55) for the C /jncomp value, results for the compressible Cf can be computed for each 
of the correlations. It turns out that results using van Driest and White/Christoph are very similar, 
so only the results of van Driest and Spalding/Chi are shown in Fig. 3 for clarity of presentation. 
Plots of Cf / C y.incomp are shown here for a variety of T w /T e ratios, as well as for adiabatic wall 
temperature. In all cases, an assumed Re x of 5 x 10 6 was used. In terms of Reg, this translates 
to approximately Reg = 14200, using the formula Reg ~ 0.0142 Re^/‘ from White [33]. For 
the adiabatic case (for which the most experimental data exist), the correlations give nearly the 
same result (black lines). But for fixed wall temperature ratios, the results can differ significantly. 
Spalding and Chi [39] claim a smaller root mean square error compared to van Driest [37] using a 
variety of experiments, but recall that relatively few experiments exist for walls with heat transfer. 
Plots of compressible Cf vs. Re x are shown in Figs. 4(a) and (b) for several specific cases. In 
Fig. 4(a), the results of the correlations agree well, whereas in Fig. 4(b) the correlations are seen to 
differ by as much as 50% or more. The main point here is that there is some uncertainty regarding 
Cf for walls with heat transfer, so it is more difficult to validate (or invalidate) turbulence models 
for these cases with confidence. In the literature, most people tend to compare with the van Driest 
correlation, but others have since developed correlations that may work better in specific cases. For 
example, Huang et al. [40] developed a method for which skin friction for strongly cooled walls falls 
below van Driest, in better agreement with data. 

The wall heat flux is often expressed [33] in terms of the Stanton number: 


St = C h = 


Qw 


(67) 


Pe C f Cj, (T,,U! 7’,,, ) 

where the heat flow at the wall q w = —k(dT/dy)\ w . The so-called Reynolds analogy is usually 
used to relate the local wall heat flux in terms of skin friction: 


St « \c f R af (68) 

where R a f is the Reynolds analogy factor. This factor generally lies in the range 0.9 < R a f < 1.3, 
and is believed to be close to unity for hypersonic flows [36] and for very cold walls [22], Thus, St 
is directly proportional to Cf, at approximately half its numerical value (with appropriate sign). 


4.2 CFD Results on a Flat Plate 

In order to test the ability of k-ui turbulence models to predict compressible boundary layer flow, 
computations were performed for flow over a flat plate in zero pressure gradient for a variety of flow 
conditions. Most of the computations were performed using the CFL3D code [15]. Note that in 
CFL3D, the turbulence models are decoupled from the mean flow equations, k is not included in 
the definition of total energy, and diffusion of k does not appear in the mean flow energy equation 
for its models tested here. Furthermore, for the current applications, Eq. (36) is used for production, 
with the (2/3 )pk term ignored. (Some computations were tried with this term included, and it was 
found to make little difference even for M = 10 cases.) The (2/3 )pk term was also neglected in 
Tij, Eq. (5). For comparison, several results were also obtained using the VULCAN code [41], in 
which the turbulence models are fully coupled to the mean flow equations and no approximations 
are made for r, ; or the turbulence production terms. 
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Table 1 . Flat plate cases computed 


Mach 

Tw/Tao 

Tyj / T au ideal 

wall type 

0.2 

1.008 

1.0 

adiabatic 

5.0 

6.0 

1.0 

adiabatic 

10.0 

21.0 

1.0 

adiabatic 

0.2 

5.0 

4.96 

hot 

2.0 

1.0 

0.556 

cold 

2.0 

2.0 

1.11 

hot 

5.0 

1.0 

0.167 

cold 

5.0 

3.0 

0.5 

cold 

5.0 

20.0 

3.33 

hot 

10.0 

1.0 

0.0476 

cold 

10.0 

10.0 

0.476 

cold 

10.0 

40.0 

1.905 

hot 


The majority of runs were performed using the Menter SST model. Full Navier-Stokes (as op- 
posed to thin-layer) was employed. For subsonic Mach numbers, the inflow boundary condition set 
total pressure and total temperature conditions (according to the particular Mach number). The pres- 
sure was extrapolated from the interior of the domain, and the remaining variables were determined 
from the extrapolated pressure and the input data, using isentropic relations. The outflow boundary 
condition set p/poo = 1 and extrapolated all other quantities from the interior of the domain. For 
supersonic Mach numbers, the inflow boundary condition set all primitive variables, and the outflow 
boundary condition extrapolated all variables from the interior of the domain. In all cases the top 
boundary, located a nondimensional distance y = 1 from the wall, used a farfield Riemann invariant 
boundary condition. The wall boundary condition enforced no slip, and set temperature either (1) 
according to a desired T^/T^, or (2) according to T w = T aw . idea i, where T awddeal is the ideal 
adiabatic wall temperature computed from T awddea i = T 00 (l + |(7 — 1 )M 2 ). This latter method 
yielded almost the same results as enforcing zero wall temperature gradient, which insured no heat 
flux at the wall. The freestream was taken to be 540 R. 

As mentioned earlier, the wall boundary conditions for turbulent quantities were those recom- 
mended by Menter [10]. Although not shown, the boundary condition on 0J wa ii was varied by a 
factor of 10 in both directions, but this change did not have an appreciable influence on the results. 
The farfield boundary conditions for turbulent quantities were determined from Tu = 0.08165% 
and Ht,oo/ Poo = (2 x 10“ ‘)Re. For the results to be shown. Re = 10 ' over the length of a plate 2 
nondimensional units long. Thus, Re per unit length was 5 x 10 6 and pt,oo/Poo = 1-0. 

The finest grid employed was 273 x 193, with 225 points on the plate and 49 points leading up 
to the plate (where symmetry was imposed). Nondimensional minimum normal spacing at the wall 
was approximately Ay = lx 10 -6 , yielding average y + = 0.2 (or less for higher Mach numbers). 
There was stream wise clustering near the plate leading edge, as shown in Fig. 5, for which only 
every fourth grid point is shown for clarity. For supersonic Mach number cases, it was necessary to 
employ a flux limiter in the computations. 

Table 1 summarizes the cases computed. For the purposes of this study, the wall temperature 
is defined as “hot” or “cold" depending on whether it is above or below the ideal adiabatic wall 
temperature. Note that the M = 2, T w /T oa = 2 case is only slightly hot, with wall temperature 
only slightly above adiabatic. 

A grid study for the M = 5 adiabatic wall case was conducted using the SST model on the 
273 x 193 (fine grid), along with two successively coarser grids for which every other grid point was 
removed in each coordinate direction (medium: 137 x 97; coarse: 69 x 49). Results are shown in 
Fig. 6. The biggest differences were near the plate leading edge, particularly for the coarse level. But 
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over most of the plate the medium and fine grids yielded very close results. For example, between 
Re x = 5 x 10 6 and the trailing edge of the plate, the Cf results on the medium and fine grids 
differed by less than 0.2%. 

In order to get an idea about the magnitude of the computed turbulent kinetic energy relative to 
the square of the velocity, profiles of local k/U 2 are shown as a function of y in the boundary layer 
at the a’ -location where Re x = 5 x 10 6 for several different cases in Fig. 7(a). The maximum level 
was only about 3%. Fig. 7(b) shows the value of another quantity sometimes used to ascertain the 
compressibility effects of turbulence (albeit most commonly for free shear layer applications) [21], 
the turbulence Mach number Mt = (y/2 ~k)/a, where a is the local speed of sound. The highest 
levels in the boundary layer at a given freestream Mach number occur for the cold-wall cases. For 
example, for M = 10, T w /T 0 0 = 1, the peak Mt reaches approximately 0.5. 

Results for the adiabatic wall cases using the finest grid are shown in comparison with van Driest 
and Spalding/Chi correlations in Fig. 8. The CFD results captured the correct trends compared to 
theory, although the tendency of numerically-induced transition in the CFD to occur further aft 
with increasing Mach number should be noted [42]. Compared to the correlations, the SST model 
slightly underpredicted turbulent skin friction for the subsonic Mach number case and overpredicted 
the correlations for the hypersonic Mach numbers. The Wilcox06 model produced similar Cf as 
SST for all three cases, although it had a greater tendency to remain laminar than SST as the Mach 
number increased. (Although not shown, increasing freestream Tu could shift the transition location 
for Wilcox06 forward.) 

Effects of a different code (VULCAN with SST model), as well as effects of including two 
different compressibility corrections, are shown in Fig. 9. The VULCAN results are seen to be 
relatively close to the CFL3D results on the same grid, yielding slightly lower Cf levels. Regarding 
the compressibility corrections, generally speaking, for adiabatic wall with freestream M < 10 there 
is only a fairly small influence on the results. Both the Wilcox and the Zeman corrections reduce the 
skin friction, with the Wilcox correction having the larger effect. Note that Wilcox [1] reported a 
larger decrease in Cf with the Zeman correction because he used different coefficients (his version 
of the Zeman model was designed for free shear flows). 

Results for two cold-wall cases are shown in Lig. 10. Again, both CFL3D and VULCAN (using 
SST) yield skin friction levels that are very close. These SST results are significantly high in com- 
parison with the correlations. It is generally well-known that simple algebraic turbulence models 
such as Baldwin-Lomax [43] can perform reasonably well for attached hypersonic boundary layer 
flows, provided that the definition of y + in the van Driest damping function uses local values for p 
and // (rather than wall values) [44], as follows: y + = \J ( pT w )y/p . As shown here, using this ver- 
sion of Baldwin-Lomax for the two cold-wall cases yields better predictions than SST, in reasonably 
good agreement with the van Driest correlation. Employing SST with the Wilcox compressibility 
correction (SST + Wilcox cc) lowers skin friction significantly. Both CFL3D and VULCAN produce 
nearly identical results. Although results still lie within the band defined by the two correlations, the 
SST + Wilcox cc results are quite a bit lower than those of Baldwin-Lomax. Results with the Zeman 
correction agree better with Baldwin-Lomax results and the van Driest correlation. As shown in 
Fig. 11, the Baldwin-Lomax model and SST + Zeman cc produce lower eddy viscosity values than 
uncorrected SST very near the wall (corresponding to y + < 0(100)). The lower levels produce a 
“less turbulent” profile, and consequently lower wall skin friction. 

Results for two hot-wall cases are shown in Fig. 12. In these cases, the results using SST gener- 
ally fell within the shaded band defined by the two correlations; for M = 0.2, T w /T 0 c = 5 results 
were closer to the van Driest correlation, and for M = 5, T w /T^ = 20 results were closer to the 
Spalding/Chi correlation. For these cases, the compressibility corrections made almost no difference 
in the results. 

Using the same plot format shown earlier in Fig. 3, CFD results using the SST model (without 
and with the Zeman compressibility correction) for various cases in terms of the ratio Cf / C'/.j ncomp 
are over-plotted alongside the correlations for Re x = 5 x 10 6 in Figs. 13(a) and (b). In all cases the 
C/, incomp used was the CFD result for M = 0.2 with adiabatic wall temperature. The SST results 
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are given by the large filled-in symbols. The trends discussed above can be clearly discerned in this 
plot. Although overall the general effects of Mach number and wall temperature on skin friction 
can be qualitatively predicted by SST with no explicit compressibility correction. Fig. 13(a) shows 
that adiabatic and cold-wall cases are consistently overpredicted quantitatively as Mach number 
is increased (up to about 30% at M = 10 for adiabatic walls, and up to as much as 40 — 100% at 
M = 10 for cold walls). On the other hand, hot-wall cases yielded results that lay in or near the band 
defined by the two correlations. Fig. 13(b) shows that SST results using the Zeman compressibility 
correction are significantly closer to the correlations. 

Results can also be plotted as a function of T w /T aw ^d ea i- Fig- 14 shows van Driest and Spald- 
ing/Chi correlations for both M = 5 and M = 10, along with current SST results. The SST results 
with no compressibility correction mostly lie above the shaded regions defined by the two corre- 
lations, except for the hot wall cases where results start to fall below the Spalding/Chi correlation. 
When the Zeman compressibility correction is employed, cold wall results are improved relative to 
the correlations, while hot wall results are essentially not changed at all. 


5 Conclusions 

The most widely-used compressibility corrections for k-u> models for high Mach number boundary 
layer flows are based on improvements intended for free shear applications. As such, these correc- 
tions are often unacceptable for boundary layer flows, and many researchers prefer to employ no 
corrections at all. Although it should be borne in mind that there is some uncertainty in the theo- 
retical correlations (especially at Mach numbers well above 5), it appears that the uncorrected k-oj 
models perform progressively worse - particularly for cold walls - as the Mach number is increased 
in the hypersonic regime. As is well-known, simple algebraic models such as Baldwin-Lomax per- 
form better compared to experiment and correlations in these circumstances. 

There is still no clarity about whether dilatation-dissipation and pressure-dilatation effects are 
important in boundary layers or not, particularly at the higher Mach numbers and for cold-wall 
hypersonic cases. Anything that reduces near-wall eddy viscosity in these situations can help obtain 
better agreement with correlations, but there is currently no strong physical argument for choosing 
one “fix” over another. Corrections designed for improving free shear flows tend to over-correct in 
the boundary layer and yield wall skin friction (and heat transfer) values that are too low. In this 
paper, it was shown that a dilatation-dissipation correction designed by Zeman specifically for use 
in boundary layer flows works reasonably well for cold wall cases. Its influence is smaller in the 
boundary layer than the popular Wilcox correction. 

The physical modeling needed to improve wall skin friction predictions in highly compressible 
boundary layer flows has yet to be formulated and accepted for k-u turbulence models. Currently, 
omitting explicit compressibility corrections works reasonably well only for lower Mach numbers 
(e.g., M < 5) or for hot-wall cases. Using the Zeman compressibility correction (formulated for 
boundary layers) improves high-Mach-number cold wall results, but it would be an insufficient 
correction for free shear flows. Better overall physics-based compressible turbulence modeling is 
clearly needed. 
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Figure 1. Contours of F c using Eq. (56). 
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Figure 2. Contours of F] j ee using (a) Eq. (62), (b) Eq. (63), and (c) Eq. (64). 
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Figure 3. Theoretical compressible wall skin friction compared to incompressible level for Re x = 
5 x 10 6 , T e = 540 R, using two different correlations. 
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Figure 4. Theoretical values of compressible wall skin friction as a function of Re x ; (a) for selected 
cases where the van Driest (solid lines) and Spalding/Chi (dashed lines) correlations agree well, (b) 
for selected cases where the correlations differ. 
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Figure 6. Wall skin friction coefficient grid study for M = 5, adiabatic wall, using SST on 3 
successive grid sizes. 




Figure 7. Profiles of local turbulent quantities in the boundary layer at Re x = 5x 10 6 , SST model: 

(a) k/U 2 , (b) M t = (v/2 k)/a. 
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Figure 8. Wall skin friction coefficients for adiabatic wall cases. 



Re x 


Figure 9. Effect of code and Wilcox/Zeman compressibility corrections on wall skin friction coeffi- 
cients for adiabatic wall cases. 
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Figure 10. Wall skin friction coefficients for two cold-wall wall cases. 



iVhnf 

Figure 11. Profiles of nondimensional //, for M = 5 and M = 10 cases, T w /T aQ = 1, at Re x 
5 x 10 6 . 
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Figure 12. Wall skin friction coefficients for two hot-wall wall cases. 


filled-in symbols = CFD 


filled-in symbols = CFD 




Figure 13. Theoretical compressible wall skin friction compared to incompressible level as a func- 
tion of Mach number for Re x = 5 x 10 6 , T e = 540 R, including (a) SST (no compressibility 
correction), and (b) SST (Zeman compressibility correction). 
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Figure 14. Compressible wall skin friction as a function of T w /T aWt id e ai for Re x = 5 x 10 6 , 
T e = 540 R, with SST results included. 
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